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ABSTRACT 

We study a maximum probability approach to reconstructing spatial maps of the Newtonian gravitational 
potential ^ from peculiar velocities of galaxies at redshifts beyond z ~ 0.1, where peculiar velocities have been 
measured from distance indicators (DI) such as the TuUy-Fisher relation. With the large statistical uncertainties 
associated with DIs (of the order ^ 20% in distance), our reconstruction method aims to recover the underlying 
true peculiar velocity field with sufficient precision to be used as a cosmological probe of gravity, by reducing 
these statistical errors with the use of two physically motivated filtering prior terms. The first constructs an 
estimate of the velocity field derived from the galaxy over-density &g and the second makes use of the matter 
linear density power spectrum Pk- Through the use of A^-body simulations we demonstrate that, with mea- 
surements with a suitably high signal-to-noise, we can successfully reconstruct the velocity and gravitational 
potential field out to z < 0.3. This will prove useful for future tests of gravity, as these relatively deep maps are 
complementary to weak lensing maps at the same redshift. 

Subject headings: methods: data analysis,numerical,statistical - galaxies: distances and redshifts - cosmology: 
large-scale structure of Universe 



1. INTRODUCTION 

1.1. Peculiar velocity fields 

The concordant gravitational instability paradigm tells us that 
the peculiar motions of galaxies arise from the effect of gravi- 
tational potentials due to inhomogeneities in the mass density. 
Therefore, these peculiar motions encapsulate information re- 
garding the growth history of large-scale structure of the uni- 
verse. Our ability to measure these peculiar velocities, either 
directly or by analysing the anisotropy of galaxy clustering 
produced by redshift-space distortions, has made their study 
an important probe of the cosmological model over the last 
two decades. 

With the advent of all-sky redshift surveys during the 1990s 
such as the 1.2 Jy IRAS (e.g. Fisher et al. 1995) and 0.6 Jy 
IRAS PSCz catalogue (Saunders 2000), the Mark III catalogue 
(Willick et al. 1996, 1997a), and the Spiral field I-band (SFI) 
catalogue (Giovanelli et al. 1998), popularity in developing 
sophisticated velocity and density field reconstruction tech- 
niques quickly grew. These methods allowed constraints on 
the bias function b and the matter density parameter, il,„ via 
the V • V = -(35 g or V = -/3V~' 5g relations, where the galaxy 
density field 5g is related to the underlying density field 6„, by 
6g = bS„i, and V~' denotes the inversion operation of V. To 
apply such relations to compute (3 one needs to have indepen- 
dent reconstructions of v and Sg (for a comprehensive review 
see Strauss and WilHck 1995). 
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Thus, central to the velocity field reconstruction method- 
ology is the acquisition of peculiar velocity data which is 
achieved through distance indicators - empirical relationships 
between two or more intrinsic properties of galaxies, that al- 
low us to estimate a redshift-independent distance d (e.g. the 
Tully-Fisher relation for spiral galaxies, or D„-cr relations for 
elliptical galaxies). In the simplest scenario, we can use this 
distance to estimate a peculiar velocity u using the following 
relation, 

u = czohs-Hod, (1) 

where c is the speed of hght, Hq is the Hubble constant and 
Zobs is the observed redshift, measured spectroscopically. 

However, the scatter in distance indicator relations intro- 
duces a typical statistical error for a given galaxy of ~'20% 
in distance at scales > 50 Mpc. Moreover, a systematic bias 
caused by the variations in large-scale structure along the line 
of sight (often referred to as inhomogeneous Malmquist bias), 
must be carefully corrected for. As a result, numerous sophis- 
ticated methods emerged such as: POTENT Bertschinger and 
Dekel (e.g. 1989); Dekel et al. (e.g. 1999), VELMOD WilHck 
et al. (e.g. 1997b); WilHck and Strauss (e.g. 1998), the unbi- 
ased minimal variance (UMV) method of Zaroubi (2002), the 
inverse Tully-Fisher (ITF) method (Nusser and Davis 1995; 
da Costa et al. 1998), and spherical harmonic decomposition 
approaches (e.g. Nusser and Davis 1994; Davis et al. 2011). 
Although all of the these approaches were extremely success- 
ful tools in probing cosmology, the large uncertainties associ- 
ated with distance indicators has thus far limited their scope 
to reconstructing spatially resolved peculiar velocity to within 
our local Universe on scales typically out to z ^ 0.015. 

Additionally, Kaiser (1987) (see also e.g. Lilje and Ef- 
stathiou 1989; McGill 1990; Hamilton 1992; Cole et al. 
1995; Dekel and Lahav 1999) demonstrated the usefulness of 
redshift-space distortions as measured directly from galaxy 
redshift maps as a more accurate probe of /3 out to deeper 
redshifts. Moreover, the advent of CMBR measurements 
e.g. made by the Wilkinson Microwave Anisotropy Probe 
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(WMAP - see e.g. Spergel et al. 2003; Larson et al. 2011) 
led to much more stringent constraints on Consequently, 
the use of reconstructed velocity fields as a competitive cos- 
mological tool has not been dominant over the last decade (for 
further discussion see also Berlind et al. 2001). 

In this paper, our goal is to apply a Bayesian maximum- 
probability approach to explore the possibility of pushing re- 
construction studies beyond these local constraints, out to a 
redshifts of z ^ 0.3. At these scales, resolved maps of the ve- 
locity field, and hence the gravitational field, could be used as 
an additional observable with which to test gravity. 

1.2. Renewed motivation 

In recent years there has been a flurry of activity to develop 
new probes of gravity on cosmological scales (see e.g. Jain 
and Khoury 2010). Several tests have emerged that involve 
cross-correlating various observables from spectroscopic sur- 
veys, galaxy imaging and the CMB. For example, it is useful 
to examine how one can constrain quantities in the FLRW 
metric in the Newtonian gauge, 

ds^ = -il-2^)dt^ + fl\f )( 1 + 2^)dx^ (2) 

where and $ are the Newtonian potentials that describe 
the temporal and spatial perturbations to the metric respec- 
tively, fl(f) is the expansion factor of the Universe and x is 
in comoving coordinates (e.g., Mukhanov et al. 1992; Ma 
and Bertschinger 1995).^ Galaxy clustering measurements 
are sensitive to the 5* component, the weak lensing shear is 
sensitive to the potentials, (^'+ $), and the CMB (via the Inte- 
grated Sachs-Wolf effect) is dependent on d{'^ + ^)/drj, where 
the derivative is taken with respect to the the comoving time 
77 (see e.g. Zhao et al. 2009, 2010; Bean and Tangmatitham 
2010). More recently, redshift distortions have been com- 
bined with other gravitational probes to further constrain dark 
energy and gravity models (e.g. Acquaviva et al. 2008; Linder 
2008; Song and Percival 2009; White et al. 2009; Guzik et al. 
2010). 

In particular, Reyes et al. (2010) combined the observations 
from galaxy-galaxy lensing (5'+$), galaxy clustering (^f) and 
the measured anisotropy in the power spectra to obtain the 
redshift distortion parameter /? as derived from Sloan Digital 
Sky Survey (SDSS) by Tegmark et al. (2006). By combining 
all these probes they were able to construct the Ec estimator 
(see Zhang et al. 2007) that provides a probe of any anomolies 
in the rate of growth of structure and/or the differences be- 
tween the potentials. Applying this estimator to the SDSS lu- 
minous red galaxies (LRGs) survey out to a redshift z = 0.32, 
they found consistency with GR on cosmological scales. 

Since peculiar velocities are sensitive to the field, it is 
of interest to investigate whether it is possible to directly in- 
fer the 4* field from velocity data at moderate redshifts of 
z ~ 0.1. This could then be compared to a reconstruction 
of the 4* + $ field from lensing on a mode-by-mode basis, al- 
lowing detailed comparison of the ^I^ and $ fields as a test of 
GR. 

The observational data is increasingly available for such a 
reconstruction (see e.g. Abate et al. 2008; Abate and Erdogdu 
2009). Recent surveys have provided considerably more com- 
prehensive catalogues of measured peculiar velocities [e.g. 
the SFIh-h- galaxy peculiar velocity survey (Masters et al. 

' In a universe without anisotropic stress, that is, where the stress-energy 
tensor is invariant under spatial rotations, or the three principal pressures are 
identical, the Einstein equations set <I> = 



2006), the 6dF galaxy survey (Jones 2009) and the 2MASS 
selected Flat Galaxy Catalog 2MFGC (Kudrya et al. 2009)]. 
There has also been ongoing work extending catalogues of 
inferred distances. These include Fundamental Plane stud- 
ies by Hyde and Bernardi (2009) which measured ^50,000 
early type galaxies in SDSS and also work by Magoulas et al. 
(2012) published 10,000 near-infrared early types in the 6dF 
Galaxy Survey out to z < 0.055. Using the Tully-Fisher rela- 
tion, studies by Mocz et al. (2012) have measured in excess 
of 25,000 SDSS galaxies at 0.045 <z< 0.085. There is also 
the emergence of another class of scaling relation to derive 
distances that relates mass and radius of Early-type galaxies 
(see e.g. Chiosi et al. 2012; Bernardi et al. 2011), which may 
have a reduced statistical error compared to current DIs. 

Several ongoing and forthcoming surveys are of great inter- 
est for potential reconstruction. The Sloan Digital Sky Survey 
(SDSS) III Baryon Oscillation Spectroscopic survey (BOSS; 
Schlegel et al. 2009) will measure redshifts for 1.5 million 
Luminous Red Galaxies over 10000 deg^ out to z < 0.6. The 
planned TAIPAN survey is intended to extend the 6df survey 
by targeting ^ 700,000 objects to give ^ 1000 objects per 
6df field. There are also the Dark Energy Survey (DES), the 
Square Kilometre Array (SKA) and the Large Synoptic Sur- 
vey Telescope (LSST) which will measure millions of red- 
shifts of galaxies on an unprecedented scale and depth. 

In this paper we revisit the use of peculiar velocity maps 
measured from distance indicators, as a cosmological probe 
for the next generation of redshift surveys to reconstruct the 
gravitational potential at redshifts > 0.1. Given the large 
statistical uncertainties inherent on the measured distances 
(^ 20%), we will apply a Bayesian approach, making use of 
two physically motivated priors to regularise and filter the re- 
construction. The first uses an estimate of the velocity field 
derived from the galaxy over-density Sg and the second makes 
use of the matter linear density power spectrum P^. Through 
the use of A^-body simulations we demonstrate that, with mea- 
surements with a suitably high signal-to-noise, we can suc- 
cessfully reconstruction the velocity and gravitational poten- 
tial field out to z < 0.3. 

The structure of this paper is as follows. In §2 we sum- 
marise the theoretical framework of density and velocity 
fields. In §3 we summarise the key distance indicator rela- 
tions and discuss some of the relevant intrinsic statistical and 
systematical sources of errors in the context of our work. In 
§4 we provide details of the A^-body simulations used to test 
the method as well as our approach to creating velocity maps. 
Details of the velocity reconstruction algorithm are presented 
in §5 together with results of applying this reconstruction to 
our simulations. In §7 we detail the method for reconstruct- 
ing the gravitational potential and discuss the results. We then 
take a brief look the prospects of applying our algorithm to 
current and future survey data in §8. Finally, in §9 we provide 
a summary and discussion. 

2. THEORY 

In this section we discuss the velocity and density field theory 
that underpins our work. 

2.1. The linear density and velocity fields 

We will consider scales on which linear theory, with zero pres- 
sure, is expected to be valid. Here the pecuHar velocity field. 
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v(x,f), and the overdensity field, (5(x,f)7 obey the dynamical 
equations for gravitational instability, Peebles (1980, 1982) 



d\ d 

-T- + -y = g, 

at a 



V • V 85 ^ 

+ 1^ = 0' 

a at 



(3) 



where g denotes the peculiar gravitational acceleration. 



g(x,0 = - 



= GpiJLa 



I over all space 

and 5* is the solution of the Poisson equation, 
V^^ = A-KGpLa^ [p(x,f)-p], 



d'-'5i-\t)^^^, (4) 



(5) 



and ^ is any modification to Poisson's equation from a non- 
GR theory of gravity (in principle this can be scale and time 
dependent, although for simplicity here we show a constant). 
The comoving position vector x is related to the physical po- 
sition vector via r(f) = a(f)x. We note that the solution satis- 
fying the second equation in Expression 3 is 



a d f g\ F(x) 
v(x,f)=-— 

47rG/i Ot \pa J a 



(6) 



where F(x) is a constant of integration that obeys V • F = 0. 
The first equation then produces 
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This is a second order differential equation which has a form 
similar to the one for 6{x, t) in the gravitational instability con- 
text, 

o2 * p) 

jd+2-—5 = 4-7:GfipS. (8) 
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dt 



Since 6 is the superposition of two modes which have a time 
dependence £>,(?) (/ G [1,2]) one therefore also expects that 



1^ cx A(0- 

pa 



(9) 



For standard FRW cosmologies, there is a growing and a de- 
caying mode. From Equation (6) one sees that the velocity 
modes are given by 



v,(x,f)^ 



gi \ Diit) 



AirGp \pa J Di(t) 
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Here, we are also assuming F(x) = implying that the velocity 
field is curl free. The growing mode will be given by 



TT r 

v+(x,f)= ^_ g+(x,f) 
HUa 



47r 



over all space 
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(12) 



HD^ ■ 70 V 2 

' 5(x,f) = p(x,/)/p— 1, where p(x,/) denotes the cosmological density 
field and p the respective mean. 



in a ACDM cosmology (Lahav et al. 1991), H{= a/a) is the 
Hubble constant, and the subscript "+" denotes the growing 
mode. Henceforward we will omit the subscript "+" as the 
only modes of we are going to consider are the growing modes 
of the different fields. 

2.2. The "traditional" approach to peculiar velocity fields 

In section 1 . 1 we briefly introduced some of the key methods 
that were developed to construct velocity fields. In this section 
we discuss these in more detail. Equation 1 1 underpins the 
"traditional" approach to draw cosmological parameters from 
a cosmological peculiar velocity field: compare v(x,f) and 
g(x,f) measured/predicted independently and compute /3. For 
example, the differential form of Equation 1 1 can be written 
as 

V-v = -/35j, (13) 

where the galaxy density field 5g is related to the underlying 
density field 5„, by 5g = bSm, and b is the putative bias function 
usually referred to as linear bias. 

Briefly, equations 11 and 13 underpin two different ap- 
proaches to the comparison between cosmological fields: the 
/) v-v and //) d-d comparisons respectively. The latter con- 
trasts the density field as measured by a full-sky redshift sur- 
vey and an estimate of the density field computed from di- 
vergence of the 3-D peculiar velocity field as derived from 
measured peculiar velocity surveys. While in the latter one 
computes the 3-D gravitational field from the galaxy distri- 
bution and the 3-D peculiar velocity field from measured pe- 
culiar velocity fields (for further detail regarding these two 
different strands see Strauss and Willick 1995). In these two 
different approaches the cosmological information is used in 
two complementary ways. In /) the 3-D gravitational field is a 
non-local field as one uses an integration over all space of the 
overdensity field see Equations 11)^ while the fields used in //) 
only incorporate local information. Our approach, described 
below, is related to //) in a Bayesian approach. 

3. MEASURING PECULIAR VELOCITIES 

Estimating an accurate peculiar velocity from current distance 
indicators (DI) remains a challenging task. In this section we 
discuss some of the key sources of statistical and systematic 
uncertainties associated with obtaining a redshift independent 
distance measurement. 

For low redshifts, the radial peculiar velocity component 
M(r) at the position r of a galaxy moving with velocity v(r) is 
given by the relation 



M(r) = [v(r) - v(r = 0)] • f = cz - //or , 



(14) 



where c is the speed of light, z is the observed redshift for the 
galaxy. Ho is Hubble constant, r = |r| is the true distance to the 
object, v(r = 0) denotes the the observer's velocity assumed to 
be at r = 0, and r represents the unit vector along the object's 
position vector r. The last equality in the above equation is 
only valid for redshifts < 0.3 (see e.g. Davis and Lineweaver 
2001), which is the regime considered in this paper If we 
are to use the relationship above to estimate u, we require not 
only the measurement of galaxy redshifts but also the infer- 
ence of a redshift-independent galaxy distance, the true value 

The overdensity field used maps the large scale structure redshift space. 
Thus the reconstruction technique also has to include a prescription to connect 
redshift distortions, which itself depends on /3 (e.g. Strauss et al. 1992; Nusser 
and Davis 1994) 



4 



R Johnston et. al. 



of which is r. The inferred distance is denoted by d and is 
obtained by means of a DI. These are empirical relationships 
between two or more intrinsic properties of galaxies, some of 
which are distance-independent (e.g. velocity dispersion cr,,, 
and Luminosity L) and some others that are distance depen- 
dent (e.g. diameter D). The observables used to characterize 
galaxies differ according to their galaxy type. Thus, distance 
indicators include the Tully-Fisher relation (Tully and Fisher 
1977) for spiral galaxies, and the Fundamental Plane (or one 
of its variations, e.g. the D„-a relation), for elliptical galaxies 
(Dressier et al. 1987; Djorgovski and Davis 1987). Both the 
Tully-Fisher and the Fundamental Plane are well established 
relations that have been used extensively since their incep- 
tion (e.g Pizagno et al. 2007; Willick et al. 1997a; WilHck and 
Strauss 1998; Willick 1999; Masters et al. 2008; Jorgensen 
et al. 1996; Graham and Colless 1997; Pahre et al. 1998; Treu 
et al. 2001; Bettoni et al. 2001; Gebhardt et al. 2003; Cappel- 
lari 2006; Springob et al. 2007; Magoulas et al. 2010; Mocz 
et al. 2012). 

The Tully-Fisher relation (Tully and Fisher 1977) is the ob- 
servation that there is a correlation between a spiral galaxy's 
absolute luminosity, L, and its circular velocity, v, as derived 
from its rotation curve, described by the power law, 

logL = A + alogy (15) 

The value of a is dependent on the fitting method adopted and 
the band in which the galaxy is observed. The Fundamental 
Plane relation (e.g. Dressier et al. 1987; Djorgovski and Davis 
1987) correlates an elliptical galaxy's surface brightness, /q, 
size, Rq and velocity dispersion, ao, such that 

logRo = alogcro + blogIe+c, (16) 

where the coefficients a, b, and c can be determined by min- 
imising the residuals from the plane by a variety of methods. 

Understanding and quantifying biases stemming from the 
above distance indicators is a complex task; here we only 
highlight the main points. A comprehensive review of the DIs 
and their idiosyncrasies can be found in the review by Strauss 
and Willick (1995). 

Peculiar velocity analysis from measured peculiar veloci- 
ties can be considered a three step process. The first step in- 
volves the calibration of the DI itself. For example, if we take 
the TF relation, calibration requires the determination of the 
zero point A, the slope a, and the scatter a. The calibration 
is sensitive to biases introduced when one uses a subset of all 
galaxies. The second step in the analysis is the inference of 
distances for sample objects. Clearly if the first step has any 
mis-calibration, there will be a bias on the inferred distance. 
Even in the absence of systematics, the uncertainties due to 
the internal scattering in the DI, crj, leads to uncertainties in 
the galaxy distances of ad/d ~ 0.2 (see e.g. Bernardi et al. 
2003; Springob et al. 2007; Magoulas et al. 2010). 

Finally, the third step is to use Equation 14 with the in- 
ferred distance d, in order to calculate the peculiar velocity. 
Note that this requires us to adopt a particular cosmology (i.e. 
the cosmological parameters assumed must be either stated or 
marginalised over). Moreover, the uncertainty in d will prop- 
agate via Equation 14 into a large uncertainty on the peculiar 
velocity for a given galaxy. From Equation 14 the error on 
the peculiar velocity is the sum in quadrature of the uncer- 
tainties in cz and Hod since these two variables are measured 
independently: 

CT« = -i/ f^cz + crjj^d ~ '^Hod ~ H^ad = HoX (0.2c/) , (17) 



assuming that typically acz will be negligible compared to 
o'Hod- Hence for instance, a galaxy at d = 150 Mpc/i"' 
has cr,, ~ 3,000 kms"' while its typical peculiar velocity is 
u ^ 300 km s~' . It is therefore necessary to average over many 
galaxies to obtain reasonable signal-to-noise. In this study we 
consider only the statistical uncertainty associated with an es- 
timated distance, convolved with a simple treatment of selec- 
tion bias, the source of which we now describe. 

With a suitable normalization, the probability of finding a 
galaxy at a real distance r can be cast as 

P(r) cx rV(r), (18) 

(e.g., Strauss and Willick 1995). Here, (j)(r) represents the 
dataset selection function in real space. Recall that we are as- 
suming that r is the real distance to an object while d is its 
estimated distance. For simplicity's sake one could suppose 
that the selection function is locally a power law of the real 
distance, say . Then, unless 2 + 9 = 0, there is a radial gra- 
dient in the probability density that leads to an asymmetrical 
random scattering of the objects when one uses the inferred 
distance rather than the real one. More specifically, for in- 
stance, if > -2 at a given distance r there are more galaxies 
perturbed from r + 6r than from r-Sr. As galaxies "carry" 
with them the velocity information at their real position, the 
average velocity field at their inferred location d will be bi- 
ased. This is the most common source of systematic error 
in the reconstructed velocity field from measured velocities 
and it is known as Malmquist bias (see e.g. Landy and Szalay 
1992; Hendry et al. 1993; Strauss and Willick 1995). Whilst 
we do not incorporate this bias fully into our analysis, we do 
include a simplified treatment which will be discussed in sec- 
tion 4.1. 

Before we describe our reconstruction method in more de- 
tail, it is convenient to discuss the simulations we will use to 
test the method. 

4. Af-BODY SIMULATIONS 

We use mock galaxy cubes extracted from the dark matter 
A^-body simulations by Warren et al. (2006) to apply our re- 
construction algorithm. These simulations adopt the standard 
concordance ACDM model, evolving initial conditions de- 
rived from the transfer functions as calculated from CMB- 
FAST (Seljak and Zaldarriaga 1996). 

In total, 16 boxes of sizes 96, 135, 192, 272, 384, 543, 768, 
1086, 1536, 2172, 2583, and 3072 /z"' Mpc were simulated, 
using the Hashed-Oct-Tree code (HOT, Warren and Salmon 
1992; Salmon and Warren 1994). For the analysis that follows 
we have chosen to apply the 384 Mpc box. 

The HOT algorithm is a particle tree-code with periodic 
boundary conditions using a Ewald scheme (Hernquist et al. 
1991). More specifically, the self-gravity of the system of 
dark matter particles is evolved according to 

2 t^t. Gm j{rj-r,\ - 

where Npart. denotes the total number of particles, egjav. the 
"Plummer softening" r, represents the position of the particle 
i-th and m, its mass. The Hubble constant H(z) for a flat Uni- 
verse at some redshift z obeys H^(z) = Hq [ri,„(l +z)-' + 51a] > 
with Hq = 100 kms"' /zMpc"' and 51a is the dimensionless 
cosmological constant. The last term on the r.h.s. expresses 
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Figure 1. Left to right: distributions of true distances (black) and our 
simulated inferred distances (red) with lognormal errors added with 
a aj = 0.3. The initial simulated cube is < zc/ 'Mpc < 384 on 
the side and we have shifted its position to 384 < zc/^ ' Mpc < 768 
(middle panel) and 768 < zc/ h^' Mpc < 1151 (right panel) to 
explore the effect of these errors at larger distances, zc denotes the 
z-component of the vector position r. 

an acceleration with respect to the position of the ;-th particle 
in the presence of the non-vanishing cosmological constant. 
The initial conditions were set up for a flat geometry universe 
with n„, = 0.3, ni, = 0.04, h = 0.7 and = 0.9. 

The galaxy velocities used in this study are found by the 
following procedure: 

• Using a Friend-of-Friends method (Frenk et al. 1988) 
we identified bound agglomerates of dark matter par- 
ticles. Only bound agglomerates with more than 400 
dark matter particles within an isodensity surface are 
considered in this analysis (Warren et al. 2006). The 
halo velocity is the average velocity of the dark matter 
particles within the halo. 

• The found halos are subsequently populated with 
"galaxies" using a model for the Halo Occupancy Dis- 
tribution of galaxies that are more luminous than a 
given minimum luminosity. More specifically, those 
halos with a mass M greater than a minimum mass 
Mmin are assigned central galaxies. These are placed 
at the halo center-of-mass and are given the mean 
halo velocity. A random number of satellite galaxies 
is then drawn from a Poisson distribution with mean 
(A^sat.) = [(M-M,™„)/M,] forM >M„i,, These satel- 
lite galaxies are assigned the positions and velocities of 
randomly selected dark matter particles within the halo. 
For more details see Kravtsov et al. (2004). In order to 
construct mock catalogues that mimic the SDSS clus- 
tering as observed by Zehavi et al. (2005) we use the 
set of values for the parameters M^in, Mi, and a given 
by Berhnd et al. (2006). 

• We construct 10 realizations of each mock catalogue by 
using different random number generator seeds when 
we a) draw a random number of satellite galaxies for 
each halo from a Poisson distribution whose mean is 
(Nsat), and b) select random dark matter halo parti- 
cles position and velocities that will be allocated to the 
satellite galaxies. The dispersion amongst the 10 real- 
izations for one mock catalogue represents the scatter 
among possible observed states for a given halo distri- 
bution and HOD model. 

4. 1 . Velocity measurement data 

Here we describe the procedure used to create our realis- 
tic noisy input velocity fields from the mock catalogues de- 
scribed above. 
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Figure 2. Example of the spatial true velocity map compared to our 
simulated observed velocity field. The top row shows l/"""^ maps at 
two different redshift slices: 0.0 < Zcos< 0.006 (left) and 0.084 < 
Zcos< 0.09 (right). The bottom panels show the corresponding ob- 
served noisy velocity maps V°^'^ where lognormal random noise of 
(jj = 0.3 has been propagated from the distances. Here Zcos denotes 
Hozc/c. 

4.1.1. Creating the observed velocities 

We test our algorithm using the simulated galaxy cube that is 
384 /;"' Mpc on the side. For each object within the mock cat- 
alogues, we have 3-D physical co-ordinates ix,y,z) and the re- 
spective 3-D velocity components (vx,Vy,v^). Since the com- 
ponent chosen to represent the line of sight is arbitrary, we 
have adopted the z-component coordinate. We will denote 
this as the real distance r of the object from the observer. As 
such the z-component corresponds to our line of sight peculiar 
velocity which is referred to as the true velocity y'""=, since 
no noise has yet been added. 

To create a catalogue of simulated observed peculiar veloc- 
ities (in the absence of selection effects) we firstly compute 
the distribution of inferred distances d given a set of true dis- 
tances r within each realisation of the mock catalogues. Under 
the assumption that errors in the magnitudes used to compute 
the inferred distances are Gaussian, the errors in the inferred 
distances are therefore lognormal. The conditional probabil- 
ity density distribution that inferred distance is between d and 
d+dd given the real distance r, is given by 

P d(lnaf) / [lnr-lndf \ 

Pr (d\r) dd = exp , (20) 

where ad is a measure of the fractional distance uncertainty 
of the Distance Indicator (Landy and Szalay 1992; Strauss 
and Willick 1995). Figure 1 shows the resulting distribu- 
tion of true (black) and observed (or inferred) distances (red) 
for the simulated cube placed at increasing distance. Here- 
after, Zc denotes the z Cartesian coordinate of the vector 
real position r. From left to right we start with the cube at 
Q <Zc /h~^ Mpc < 384 (corresponding to maximum redshift of 
[HqZc/c) ^ 0.018; hereafter Zcos denotes Hqzc/c) and observe 
that the effect of the errors only becomes significant beyond 
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Figure 3. The correlation r between the true velocity field V'™"^ and 
observed velocity field V"'"' for ad = 0.3 at increasing redshift slices. 
The fitted slope is shown as a red line and given by the quantity a. 

Zc//!~'Mpc > 200. We then shift this cube in distance by its 
maximum dimension, 384 /;"' Mpc to 384 < zc//i~'Mpc < 
768 (middle panel) and 768 < zc/h'^Mpc < 1151 (right 
panel). We can clearly see that beyond 300 /;"' Mpc the orig- 
inal zc distribution has been substantially smeared. 

4.1.2. Creating pixelized velocity fields 

We now wish to average the data on a particular scale. This 
is achieved by pixelising the cube equally in x, y and redshift 
ZcosWith a transverse pixel scale of (24 Mpc)^. Thus we 
have 16x 16x 16 = 4096 pixels in total. This scale is suf- 
ficiently large to smooth out non linearities in the velocities. 
For this analysis we consider the simulated galaxy cube at two 
redshift ranges: 0.08 < z < 0.18 and 0.18 <z< 0.27. Their 
respective number densities are summarised in Table 1 . The 
first column in the table shows the redshift range of each sim- 
ulated cube used in our analysis. The second column shows 
the pixel scale in square degrees. For each simulation the pixel 
scale is computed based on the midpoint of the redshift range 
as indicated. The third column shows the total number of ob- 
jects in the catalogue and the final column is the number of ob- 
jects per square degree. For the 0.08 < z < 0. 18 cube we have 
a total of 303, 575 objects which equates to an approximate 
number density of ^ 140 galaxies per square degree. Shifting 
this cube to 0.18 < z < 0.27 corresponding to number density 
of ^ 341 galaxies per square degree. For comparison we in- 
clude the estimated number density from the SDSS Legacy as 
detailed in SDSS-III Collaboration: Christopher P. Ahn et al. 
(2012) at a redshift range of 0.05 ^ z ^ 0.15 comparable to 
the 0.08 < z < 0.18 simulation. These details will be used in 
the analysis in section 8. 

We bin the data in redshift space and for each slice in Zcos, 
we can evaluate the average velocity Vijk for the (/, j)-th pixel 
in the ^-th redshift sUce 



J^uiW {rijk-r,\L) 

j:w{rtjt-r,\L) 

1=1 



(21) 



Table 1 

Number density of objects within the simulated galaxy cube. For 
comparison we include details of the SDSS Legacy survey for a 
similar redshift range as the 0.08 < z < 0.18 simulation. 



Redshift range 


pixel scale 






z 


(sq. degree) 




(per sq. degree ) 


Our simulation: 








0.08 <z< 0.18 


8.47 (at z = 0.13) 


303,575 


~ 140 


0.18 <z< 0.27 


3.47 (at z = 0.225) 


303,093 


~ 341 


SDSS Legacy: 








0.05 <z< 0.15 


8032 (total area) 


~ 536,000 


~66 



144,660t 



' Approximately scaled to the 0.08 < z < 0.18 simulation 

where A^g^i is the total number of galaxies in the catalogue and 
M/ is the /-th velocity of the galaxy within the catalogue, and 
the kernel W is the 3-D top hat given by 

W {rijk - r, \L) = Wixij,-x, \L) x Wiyijt -y, \L) 

xW(zijk-zi\L), (22) 

where 

{1, |jt;| < L/2 or x = L/2 
(23) 
0, otherwise. 

In these equations r,jA is the center of the (/, j, k)-th pixel 
while r/ is the vector position of the Z-th particle. In the cur- 
rent analysis we consider velocity field reconstructions at each 
slice along the zc-axis for the true and observed velocity fields 
denoted respectively as V'""^ and y*^ Figure 2 shows an ex- 
ample of the resulting velocity maps for two redshift slices 
at 0.0 < Zcos< 0.006 and 0.084 < Zcos< 0.09 to demonstrate 
the effect of the noise term in the V°^'^ maps. The top panels 
represent the true velocity field V'""^ at these redshifts, while 
the bottom panels show the corresponding observed velocity 
fields y°''^ The labels 'X' and 'Y' on each panel represent 
the x-y pixel positions in real space. The z-axis shows the 
magnitude of the velocities in km s"' . 

As expected, for the nearby redshift slice on the left hand 
side of the figure, there is little difference between the V'™'^ 
and V°^'^ maps. However, as we move to the highest red- 
shift slice (right-hand panels) in the simulated cube at < 
Zc//!~' Mpc < 384, the original velocity map signal is com- 
pletely noise dominated at the ctj = 0.3 level. This is clearly 
shown in Figure 3 where we examine the correlation between 
the [yt™e yobsj distributions. In each panel we estimate the 
correlation between the two quantities by applying the sample 
coefficient of correlation r(y'™, V°^'^) and a standard least- 
squares approach to estimate the slope a which is also indi- 
cated by the red line, with a one-to-one relation shown as a 
grey dotted line for comparison. For the < Zcos < 0.006 slice 
we see a very tight correspondence between the [ytrue yobsj^ 
with r = 0.981;!^ oQ5. As we move to higher redshift ranges 
we can see how the noise begins to dominate; at 0.174 < 
Zcos< 0.180 there is very little residual correlation between 
[y"™,y°''^] with r = 0.174! 



-1+0.116 
-0.121- 
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Table 2 

Main input parameter values for reconstruction algorithm. 



the likelihood: 



Parameter 


yrsx 




N, 


15 


15 




20000 


20000 




1000 kms"' 


1000 kms-i 


'^m 


0.1 


1000 s-2 
0.1 



In this section we describe our velocity reconstruction al- 
gorithm. The task of a reconstruction algorithm is to infer a 
good estimate of the underlying true quantity T. The inferred 
quantity, T, can be found by minimising misfit statistics tak- 
ing into account the noise of the observations, A^, and the data 
D. The particular statistics that one applies depends on the na- 
ture of the noise. Having chosen our goodness-of-fit or mis- 
fit statistics a common prescription is to parameterize T and 
try to extract such parameters from the observed data, D, via 
B ayes' theorem 



Pr(f 



Pr(D|f)xPr(f) 
Pr(D) " 



(24) 



Pr(D|7') is the likelihood while Pr(7') denotes the prior, which 
contains any previous information or prejudice we have re- 
garding the value of T. In our algorithm, the relevant data are 
the noisy measurements of the radial velocity field V°^'^ and 
the quantity we seek to infer is the reconstructed radial veloc- 
ity field y'"'^'^-. This is an application of maximum probability 
reconstruction to velocity distance measures; much work has 
been done in applying this methodology in lensing e.g. Seitz 
et al. (1998); Bridle et al. (2000); Marshall et al. (2002) and 
to large-scale structure e.g. Jasche et al. (2010); Kitaura et al. 
(2009). In our approach we introduce a set of physically moti- 
vated priors crucial to filtering the velocity reconstruction and 
aiding convergence toward the true underlying velocity field. 
In this context with a particular choice of priors, the previous 
equation can be cast as 



log[Pr(y'' 



|y°bs)]=log[Pr(y°bs|y: 



)]+ 

+ log[Pr(Pv(k)|P*(k))] + 
+ log[Pr(y|(5p]-const., 



(25) 



where Pv(k) is the measured velocity power spectrum and 
Pyik) is its theoretical conterpart. The quantity 6^ is an ex- 
pected velocity field derived from the galaxy overdensity. In 
essence, our choice of priors assumes linear theory, where the 
Sg prior is constructed directly from the galaxy distribution 
via the Poisson equation. The two priors compliment each 
other with the power spectrum regulating the likelihood on 
various spatial scales, and the (5j field filtering on a pixel by 
pixel basis. The construction of the priors will be be discussed 
in detail in the following subsections. 

The second and third terms of the previous equation are the 
logarithm of the prior while the first one is the logarithm of 



-2iog[Pr(y*^|y--)] = x,'= 5] 



pix.k /yrec _ T/obs\2 



(26) 



where A^pix.jt is the number of pixels in the A:-th zc slice and 
'^lijk is given by 



!:'(«, -y*^)V2 (r,,,-r,|L) 
/=i ^ ^ 



(27) 



(r,v,-r,|L) 

/=i 



The total number of objects in a given pixel (i,j,k), Nijk, is 
given by 



/=1 



-r,|L)=^W^^i 

/=i 



rijk-ri\L) 



(28) 



V°''' is the average velocity in each pixel as defined in Equa- 
tion 21. 

We note at this point that in its current form. Equation 26 
is going to a biased estimator; the x^distribution is correct 
for data with Gaussian errors; however, as shown in Expres- 
sion 20 our inferred distances d are drawn from a skewed, 
lognormal distribution which is subsequently propagated into 
the observed velocity field V°^^. We can, however, introduce 
a corrective term into the likelihood to account for the mean 
offset caused by this error distribution. Recall that for each 
object / we have defined an inferred distance c/, and a true 
distance r, which were then converted into effective redshifts 
Zcos.i and Zobs.i as Hodt/c and Horj/c, respectively. Thus for a 
given slice in redshift we define a correction A'^"" as 



(29) 



Figure 4 shows how A™" for the 0.08 < z < 0.18 simulated 
cube varies with redshift. We compute A*^""" by generating 100 
different catalogues of observed distances and computing the 
(Zobs~Zcos) Izok ^iid dispersion on this for each realisation. The 
red line shows the least squares fit to the realisations. With the 
inclusion of the corrective term A'^°" the log likelihood term 
is given by: 



Xk 



2_ ^ ^^Uk ^7^+^001x1)^ 



(30) 



'■j=i 



vJjk 



Having obtained a form for the likelihood, we now turn to the 
prior terms. 

5.1. The 6y Field Prior 

In this subsection we discuss the first prior term used in 
the reconstruction. We construct an additional velocity field 
c)^'-; derived from the observed over-density (5^ of the galaxy 
distribution. We first compute the over-density 5g along the 
line of sight 



N-N, 



(31) 



where is the number of objects in each pixel at a given z 
slice and A^, is the mean number of objects within each pixel 
at a given z slice. From this starting point we can construct an 
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Figure 4. Correction applied to the velocity x term, to remove offset 
of observed distances from true distances. 

estimate of the gravitational potential field ^ via the Poisson 
equation (assuming a hnear bias b = 1). In Fourier space the 
gravitational potential is given by 



(32) 



where Vl„, is the dimensionless matter-density parameter, Hq 
is the Hubble constant and k is the wave-vector The scale 
factor a{t) is computed at the midpoint of each redshift slice 
given by 

a{t)=-^. (33) 

The estimated velocity field d"^, is obtained by taking the 
partial derivative of along the the line of sight. In A:-space 
this is computed by 



(34) 



Here is the z-component of k. Therefore this (5,, prior term 
takes the form 

;=1 



5,.= 



(35) 



The error ct^, is a parameter that determines how much this 
prior term will dominate in the reconstruction. The value we 
adopt is given in Table 2; the implications of this will be dis- 
cussed in more detail in section 6.2. Writing the term in this 
form assumes that off-diagonal terms in the covariance matrix 
are small. Further tests have shown this to be the case. 

As a quick consistency check we can compare the 9^',- field 
with the true velocity field V'™*^ to see how well they are corre- 
lated. Figure 5 shows the [V'™*^, d'i> ] for a sample of redshift 
slices < z < 0.006 (top left), 0.017 < z < 0.023 (top right), 
0.129 < z < 0.135 (bottom left) and 0.174 < z < 0.180 (bot- 
tom right). In each case we find a high correlation between 
the distributions with marginal scatter We observe slight off- 
sets in the amplitudes of the velocities indicated by the slope 
of the red line. 
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Figure 5. Distribution of the velocity field V""" and the 93'; field 
derived for selected redshift slices. The quantity r is the correlation 
coefficient. The quantity a is the slope from a least-squares fit with 
the best-fit shown as a solid red line. The one-to-one relationship is 
shown as a dotted grey line. 

Figure 6 shows a spatial representation of the steps to cre- 
ating the field for the < z < 0.006 slice. The first (left) 
panel shows the number of galaxies per pixel (A^gai); the sec- 
ond panel is the field as computed from Equation 32; the 
third panel is the 9^1'; velocity field as computed from Equa- 
tion 34; and finally the true velocity field V'™*^ is shown in the 
fourth panel. The features of V and 9^1', are very similar. 

5.2. The Peculiar Velocity Power Spectrum Prior 
The second prior term makes use of the theoretical velocity 
power spectrum to regularise the V'^'^'^ field. Using this infor- 
mation in reconstructions is not a new concept; Kaiser and 
Stebbins (1991) demonstrated the usefulness of adopting a 
power spectrum prior for regularising non-parametric fits in 
their work on density field reconstruction from peculiar ve- 
locities. The (5,, prior from the previous section compares the 
velocity fields on a pixel by pixel basis. As we shall see, the 
introduction of the P/; prior complements this by regulating 
the reconstruction on a variety of scales. 

In order to use this prior, we need to compute the power 
spectrum of the true velocity field V'™*^ and the trial recon- 
structed velocity field V"^'^. We compute the power spectrum 
over the entire simulated cube in the following way. For a 
given simulated cube we Fourier transform the field, and ex- 
amine spherical shells in wavenumber For each shell we 
compute a velocity power spectrum as 



N,(k) 



(36) 



where Na{k) is the number of pixels in each shell for each 
^-mode. A^pix runs through all the pixels in the cube. The 
associated error is given by 



[PAk)^-PAkf]'/^ 



(37) 



VNa(k)-l 

We use the true velocity power spectrum P'J^'^ik) as a proxy 
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Figure 6. Comparing the S^*; and V'""^ velocity spatial maps of the redshift slice < z < 0.006 from Figure 5. First panel shows the total 
number of galaxies within each pixel for that slice. Second panel: from this we compute the field from Equation 32. The third panel shows 
the resulting gradient of the "^^"^ field, which we can then compare to the true velocity field V'™"^ shown in the fourth panel, computed from the 
z-component velocities of the simulation. 



for a theoretical power spectrum. We compare this theory 
velocity power spectrum with each trial reconstruction 

power spectrum so that the Pk prior takes the form. 



E 



true 1 2 



[p,iky<"'-p,{k m 



(38) 



where k^ot is the total number of wavenumbers and 

f7];(,) = flP„W™. (39) 

The parameter a in Equation 39 allows us to modify the strict- 
ness of this prior term in the reconstruction. We have found 
that the off-diagonal terms are small and therefore can be ne- 
glected in the above prior term. 

5.3. Implementing the algorithm 

We will now describe the implementation of the algorithm. 
Figure 7 shows a flow diagram which summarises the follow- 
ing details: 

Step 1: First we create A', initial trial velocity fields (we 
choose N, = 15''), which we caU V™'}"' (FFT denotes fast 
fourier transform, i.e. we are working in Fourier space). We 
do this by simply taking the true velocity field V^'fj and add 
gaussian random noise with a dispersion, rj, that controls by 
how much we perturb the field. Variables x*' and are ran- 
dom numbers drawn from a gaussian distribution with mean 
zero and variance of unity, where is the imaginary part. 
Although, in practice, one would begin instead with a per- 
turbed field, tests have shown the final result to be fairly 
robust to the initial guess. Clearly, the larger the value of r/ 
the greater level of noise added and thus the worse our initial 
guess becomes. 

Step 2: Once we have established our initial trial velocity 
fields we begin the main loop of the algorithm: for each itera- 
tion of the reconstruction we perturb each n = 2,Ni trial recon- 
struction by adding to each mode an amount of uniform ran- 
dom noise (x" and y" are random numbers drawn from a uni- 
form distribution between [-0.5,0.5]), the magnitude of which 
is governed on a mode by mode basis by the true power spec- 
trum Pl'"'^(k). The « = 1 field is not perturbed (see step 5). The 
parameter A^u is a cooling term that regulates the amount of 
pixel variation allowed for each iteration and thus aids con- 
vergence. This is discussed in more detail in Step 5. 

Step 3: For each «,/, trial reconstruction V^pj{n) we com- 
pute its corresponding power spectrum P^'^^in). We also trans- 

' The choice of Ni is arbitrary. Tlie value we selected was motivated by 
code optimisation and to allow enough variation of trial fields 



form each Vppj(n) to real space for evaluation in the log like- 
lihood. 

Step 4: Compute the log probabilities for each n trial recon- 
struction. 

Step 5: Determine the minimum log probability from the fol- 
lowing equation for the A', trials 



X («) = > , (40) 



1=1 

[V""'(n)-d'^(n)]^ 

i=l 

[Ff'=(n)-P'™(n)]2 



-E 



'p(k) 



where n is the n''' trial reconstruction. If the trial with mini- 
mum log probability corresponds to « = 1, we resample from 
the parent sample in the next iteration of the code and allow 
the pixel variation parameter Ap,j. to reduce for the new set of 
trial fields. However, if the minimum log probability corre- 
sponds to « y 1, we set the corresponding trial V'^^'^ to be the 
new template and feed this back in to Step 2. Steps 2 to 5 are 
repeated until convergence has been achieved. 

6. RESULTS 

We have applied our reconstruction technique to a simu- 
lated galaxy cube at two different redshift ranges: 0.09 < z < 
0.18 and 0.18 < z < 0.27. It should be noted that the we have 
used the same set of galaxies for both redshift ranges; this 
provides a useful test of how the reconstruction changes as we 
move to higher redshifts. We have performed a Monte Carlo 
(MC) analysis to explore the error distribution of the veloc- 
ity reconstruction. For each MC we maintain the same initial 
input observed velocity map, but we allow a different ran- 
dom seed in Steps 1 and 2 of Figure 7 where we perturb both 
the initial trial velocity fields and the trial reconstructed fields 
respectively. We have performed 1000 MC reconstructions, 
which is adequate for our error analysis below. For each cube 
we have applied identical input parameter values, summarised 
in Table 2. The a regularisation parameters of the priors are 
selected carefully so that the prior terms do not dominate the 
reconstruction (see discussion in section 6.2). 

Figure 8 shows final results of selected slices from the two 
simulated cubes. In the upper panel sets we show the result 
of the reconstructed velocity field for the last redshift slice 
(0.174 < z < 0.180) in the 0.09 < z < 0.180 cube. On the 
left-hand side we show the true velocity map V'""^ in the top 
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set up Nt initial trial velocity fields 

yjP^y (n) = Vp'^^ + Vp^^^jx^ + iyg)77^ 



Begin loop at n=2 to retain previous field 



Step 2 



perturb initial fields for trial reconstructions 

VpFTin) = Vptkn > 1) + [x^ + zy'^)/P(fc)^A 



Step 3 



For each trial reconstruction compute its power spectrum P^^^{n) 

( Perform inverse FFT of each trial y^'py(n) reconstruction | 
i for the X calculation j 



Step 4 



End of loop 



i=l 



i=l 



+ 



'p(k) 



Step 5 ( Find min[;t:^(n)] j< 



if min[;^^(n)] is(n = 1) 


, resample from (n=l) and reduce ^pix ^ 




if min[;^^(n)] is (n > 1) 


,set F*"^i(n = l,Nt) = F(niin[x^]) ] 



Perform FFT of new trial velocity field and 
feed back into Step 2 for next iteration J 



Figure 7. Flow chart summarising the velocity reconstruction algorithm. See main text for details. 



row, followed by the reconstruction V'^'^'^ in the middle row. 
The reconstructed map is the mean of the 1000 MC sam- 
ples. The two panels in the bottom row show the correlation 
of r[V'™^y''='=] (bottom-left) and r[y'™^y°'^^] (bottom-right). 
The slope a of both distributions is shown and also indicated 
by the red line. 

These plots reveal promising reconstructions, comparing 
well to the underlying true velocity map with a high cor- 
relation coefficient value of r = 0.848;!;o and a slope of 
fl = 1 .000 ± 0.04. This is promising when we compare this to 
the input observed velocity map, which showed only a small 
residual correlation between the observed and true velocity 
fields, r = 0.123;!;q [23. One could suppose that such a good 
reconstruction from noisy data is purely the result of a dom- 
inating (5,, prior, but as we will demonstrate in the following 
section that this is not the case. 

If we turn now to the corresponding right-hand panel set 
for the 0.174 < z < 0.180 slice, we show the residual map 



denoted as e, which is defined as e = (y""'= - V'^^^). The mid- 
dle row of this set shows the standard error cr for each pixel 
and is computed as follows. For each pixel located at (i,j) in 
each redshift slice we compute the mean V'^'^'^ of the result- 
ing binned velocity distribution of the 1000 final reconstruc- 
tions. The right-hand panel of Figure 9 shows the distribution 
of MC reconstructed velocities for a typical pixel. The left- 
hand panel of Figure 9 shows an example of the distributions 
of correlation coefficient values we have determined for the 
MCs. We estimate the standard error on the velocity in the 
pixel as 



o-i.i = 



yrec] 



(41) 



The error cr map for the reconstruction shown in Figure 8 rep- 
resents a typical range of a that is ~ [85,100] kms"'. This 



Velocity reconstruction of gravity 



11 





Figure 8. Final reconstruction of the 0.09 < z < 0.18 (upper panel sets) and 0.18 < z < 0.27 (lower panel sets) simulated galaxy cubes. For the 
left-hand panel sets, the top row shows the true velocity field and the middle row is the final reconstruction from the mean of 1000 MCs. The 
bottom row left panel is the distribution of and the bottom row right panel is the observed velocity distribution V°*"' compared to the 

true velocity distribution V'™"^. In both of these panels we show the correlation r between the variables and the fitted slope a (also indicated by 
the red line). The right hand panel sets show the residual map (top row) e =v""'^-V"^'^, the MC error map (middle row) and the covariance of the 
reconstructed velocity map (bottom row). 
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Figure 9. Left panel: An example of the distribution of the correla- 
tion coefficient r for the 1000 MC velocity reconstructions (Vmc) and 
the true velocity field (V'™") at a redshift slice 0.264 < z < 270. Right 
panel: Velocity distribution for 1000 MC reconstructions in a single 
pixel. In both panels the results were derived from the second sim- 
ulated cube which was shifted to the redshift range 0.18 < z < 0.27 
with a (7rf=0.3 error added to the observed distances. 

is clearly considerably lower than the (T5,,=1000 kms"' value 
set as the regularisation value for the 6v prior term and is thus 
indicative that this prior is not dominating the reconstruction 
process. 

Finally, the bottom row shows the covariance of the recon- 
structed velocity map which, for a given slice in redshift. We 
find very little covariance between the pixels. This is not sur- 
prising since the large pixel scales (24 /!"'Mpc)- over the area 
we are reconstructing in each redshift slice should show little 
correlation. 

In the bottom panel set of Figure 8 we display the re- 
constructed velocity field for the 0.264 < z < 0.270 slice in 
the 0.180 < z < 0.270 cube. Once again, even out to this 
high redshift, we are still able to recover a reasonable ve- 
locity map with r[y"™,y'^'='=] = 0.853!;j;j]^° and a slope of 
a = 1.294 ±0.052. Visually, in the 3D-maps of and V''"' 
we can see a suppression of the reconstructed velocities com- 
pared to that of the true velocity field, however, it is clear that 
the overall structure has been recovered well. 

In Figure 10 we show the overall trend of the correlation 
between V''^'^'^ and V^'"' for each slice in the simulated cubes. 
In the bottom panels of the figure we show correlation coeffi- 
cient values r between [ytrue^yobs] (red points) and [V^^y 
(black points) for the simulated cube at 0.08 < z < 0.18 (bot- 
tom left) and 0.18 < z < 0.27 (bottom right). In the upper 
panels we show the corresponding correlation coefficient val- 
ues between the residuals e and: V'™'' (black points) and W'"^ 
(blue points). It should be noted that the line through connect- 
ing the points is not fitted, but only there to aid the eye. 

We observe that the overall correlation with the input ob- 
served velocity map V''^'^ across the two cubes is predictably 
low, with a maximum r = 0.228;!;Q j|9 and a minimum r = 
-0.008 1;!;"' [22- This is due to the addition of our large <Jd = 0.3 
error on the distances at these relatively high redshifts. 

In contrast, the black points on the lower panels of Fig- 
ure 10 demonstrate that we achieve a good reconstruction 
within each redshift slice. We note that the excessive dip in 
correlation in the 0.18 < z < 0.27 cube within three consec- 
utive slices at: 0.197 < z < 0.202, 0.202 < z < 0.208 and 
0.208 < z < 0.214 is due to two outlying pixels with ex- 
tremely high velocities of the order ^6,000 kms"' which the 
algorithm has failed to constrain effectively. 

We will now examine aspects of the reconstruction in more 
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redshift (z) 

Figure 10. Correlation coefficient values for the final reconstruction 
of the peculiar velocity field. The bottom panels show the corre- 
lation between r[l/°''\K'™'^] (red points) for each redshift slice in 
our simulated galaxy cubes: 0.09 < z < 0.18 on the left panel and 

0. 18 < z < 0.27 on the right panel. The black points show the mean 
correlation value for the 1000 MC samples and the true velocity field 

1. e. r[ V"^'^,V'™°]. The top panels show the respective residual corre- 
lations between the (black points) and V"^'^ (blue points), where 
e is the residual, e = V""" - V""". 

detail. In particular, we look at how the velocity power spec- 
trum of the reconstruction evolves throughout the reconstruc- 
tion process and compare it to the power spectrum of the 9^. 
and y'""^ fields. We then examine how the priors affect the 
outcome of the reconstruction when individually applied. Fi- 
nally, we examine how observational effects impact the con- 
struction of a field and consequently how the (5,, prior will 
bias the final reconstruction. 

6.1. The residual power spectrum 

A useful statistic with which to evaluate the reconstruction 
is the power spectrum itself. In Figure 1 1 we have plotted 
the power spectrum as measured from the observed velocity 
field V"^"'' shown in green, the true velocity field shown in red, 
the 95' J. field shown in blue, and finally the residual power 
= Pi(95'z- y"^'^) after the first n = 1 iteration of the code 
(dashed black line) and the final n = 20,000 iteration. 

The difference between the P*'(/t) and Pj™(/t) is self ev- 
ident, with the observed power being flat across most scales 
due to the dominant measurement noise, f '™'^(fc) and P^^{k) 
are rather similar, showing low power on large scales. The 
residual power ep^ reveals important information on how the 
reconstruction evolves; if we look at ep, after the first iteration 
of the reconstruction then we can see that already the shape of 
its power spectrum is similar to that of Pl"'^(k) and P^^{k) but 
with a comparative excess on all scales. By the final iteration 
of the code « = 20,000, ep, (black solid line) now has similar 
power to P^'^(k) on small scales k < 0.05, but less on scales 
larger than this. At k> 0.05 we observe a similar trend when 
compared to P^''''{k). This is a promising result as it shows 
that the residual power is smaller than the true fluctuations on 
these scales. 
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Figure 11. The error power spectrum compared for the various veloc- 
ity fields. The red line shows the power spectrum of the true velocity 
field V'""^, the blue line is the power spectrum on the d^- field (com- 
puted from the galaxy over density) and the green line is the power 
spectrum of the observed noisy input velocity map V°''\ The dashed 
black line shows the residual power of (d^z - V"^"^), where V"^'^ is 
the reconstructed velocity field for the first n= I iteration of the al- 
gorithm. The solid black line is then the residual power for the final 
iteration of the algorithm at n = 20, 000. 

6.2. The impact of the priors 

All the results presented so far have used the power spec- 
trum P{k) prior term and the 5y prior term in combination. In 
this section we examine how each prior is affecting the re- 
construction process individually. In Figure 12 we show four 
scenarios of possible reconstruction processes for one of the 
MC realisations at 0.163 < z < 0.169. In all cases we show 
the [y'='^,y'™'=] distributions with their respective correlation 
coefficient values r and the slope of their respective slopes 
a. Firstly, we test the reconstruction method where no pri- 
ors are activated as shown in the first panel (left hand side) 
of the figure. This essentially serves as a null test, where we 
are essentially fitting to the observed velocity data which is 
a convolution of the true velocity information and the log- 
normal noise derived from the simulated observed distances. 
As expected, the plot shows a very noisy reconstruction with 
a correlation of r = 0.079;!;o [j' and a very broad distribution 
with velocities in excess of ^3000 kms"' at the tails. 

In the second panel we show the reconstruction where only 
the prior has been activated (with the same value for crj^.^^^ 
as used previously, given in Table 2). What is immediately 
clear is a much improved tighter distribution between V'™'' 
and V'""' with a correlation of r = 0.420:l:„ 'iog. We also find an 
improved slope a of the correlation, a = 0.38 ± 0.05. 

The third panel of Figure 1 2 shows the effect of the 5y prior 
which compares the V"^'^ field with 9^. field derived from 
the galaxy over density. We adopt the same ct^, =1000 kms"' 
value as in our previous analysis. The motivation for choosing 
such a conservative error on the prior is to limit its effect on 
regularising the reconstruction and avoid us simply fitting the 
field. As we can see, the final reconstruction with the 



5y prior alone shows an improved but still low correlation of 
r = 0.39 1!|]:°(5^ and a slope of a = 0.2 ± 0.03. It is clear that our 
choice of as^. allows for a mild filtering of the reconstruction, 
but not to the extent of dominating it, which would lead to 
r ~ 0.8 as in Figure 5. 

Finally, by combining the two priors as shown in the 
fourth panel, we find a much improved reconstruction with 
r = 0.661!jj;j]^^ and a = 0.63 ±0.04. This demonstrates the 
importance of the combined effect of these priors for a suc- 
cessful reconstruction of the velocity field. 

6.3. Observational effects 

Applying our method to real survey data will, of course, be 
subject to various selection effects and biases. Understanding 
how this might impact the 6y prior, and consequently our ve- 
locity reconstruction, is necessary if the method is to be suc- 
cessfully applied to future data-sets. In the final part of this 
section we test how much the 5,, prior will bias the reconstruc- 
tion if in some region mass does not follow galaxy counts, i.e. 
if there is a baryonic feature present which is not represented 
in the true gravitational potential field. We test this by arti- 
ficially embedding into the gravitational potential field 'fake' 
over- and under-dense regions which are not in the true veloc- 
ity field. We calculate the resulting d'i>, field and hence the 
(5,, prior term. The purpose of this test is to observe if such 
effects will substantially contaminate the final reconstruction 
of the velocity field. 

6.3.1 . Embedding an overdensity 

In Figure 1 3 we consider the scenario where a large mass 
would be infeiTed to be present from the galaxy clustering, 
modelled by a gaussian embedded into the field over sev- 
eral redshift slices. In the figure there are two panel sets that 
show two of the redshift slices where this false mass is most 
prominent. In the left-hand panel set (0.135 < z < 0.141) 
there are two columns of results: the left-hand column is the 
reconstruction process for that slice under ordinary initial in- 
put conditions, where the 5*°^' field has not been artificially 
altered. The right-hand column shows the reconstruction re- 
sults where a fake mass has been embedded in to ^1*°^'. To 
distinguish the two we will refer to the altered prior infor- 
mation with a "hat" above the observable e.g. Vl/°^'. The top 
row shows the two ^E*^^' maps where a noticeable peak in ^P"^' 
map can be seen centred at [11,12] on the [X,Y] pixel plane. 
This is then propagated into the field where we now ob- 
serve negative excess in the velocity potential extending to 
^ 600 kms"' centred at the same pixel position (second row). 
The difference between d'i'- and d'i>, is perhaps more clearly 
shown in the third row where we plot the distributions of both 
[9*,,y'^"'=] and [9*,,y'™]. In the right-hand panel of this 
row we now see a skewed distribution of velocity pixels to- 
ward the upper right part of the plot compared to the d'i^, on 
the left where it remains a relatively tight distribution. 

Our new d'i'r field provides a modified S^, prior term of the 
form 

(42) 

We use our usual reconstruction methodology with this prior 
and the standard Pk power spectrum (our regularisation and in- 
put parameters remain unchanged from Table 2). The bottom 
row of Figure 13 shows the distribution of [y'^'='=,y'™'=] (left) 



14 



R Johnston et. al. 



3000 




|-a= 0.01 ±0.01, ■ 
0.163<zg0.ie9 



5 -3000 - 



T — r 



1 — r 



I U.U f »_o ,230 



T 



no priors 
I . . I . 



_L 



-3000 3000 



-1000 - 



I '-'■^<'"_0.1084 



r=0 riQI +01991 

I J- -0.1091 




r = 661 +""^3' 

I "J "" i -0.0752 



[<5„ + P^l priors" 
• I ... I 



-1000 1000-1000 1000-1000 1000 



Vrec (kms->) 



rec 

Figure 12. The effect of the different prior terms on the reconstruction. The left panel shows the final reconstruction of the velocity field where 
no prior terms have been activated. The second panel shows the final reconstruction when only the power spectrum prior term is applied with 
a (7)i(j.)=0. 1P'™'^(A:). The third shows the reconstruction with the baryon 5,- prior term only with a (7a,, = 1000.0 kms"'. Finally, the fourth shows 
the combination of the two prior terms. 



and [y'"'='=,y™'=] (right). Comparing the two distributions, we 
can clearly see that there is very minimal residual bias in the 
final reconstruction of the velocity field despite the prominent 
artificial velocities present in the 95', field. 

The left-hand panel set of Figure 1 3 shows the same recon- 
struction for the adjacent redshift bin at 0.141 < z < 0.146. 
Once again, we find the final reconstruction to be marginally 
adversely affected by the presence of the embedded mass. 
These results are indicating that the presence of the prior 
term and the velocity likelihood term are suppressing the ad- 
verse observational affects in the 9*1'; field that should not be 
present in the underlying true velocity field. Importantly, this 
result demonstrates that the (5,. prior term does not dominate 
the reconstruction process. 

6.3.2. Embedding an under density 

We now embed an additional under-dense region of a simi- 
lar magnitude to that of the over dense region in the previous 
section. This represents the case where the galaxy distribution 
under-represents the gravitational potential in a region. We 
keep the same artificial mass from the previous test and now 
add the additional under-density over different redshift slices. 
In Figure 14 we have again selected redshift bins where the 
feature is most prominent i.e. 0.163 < z < 0.169 (left panel 
set) and 0.169<z<0.174 (right panel set). Compared to the 
case where no under dense region has been (third row down 
in Figure 14 we observe an excess region of velocities toward 
the right hand region of the relevant plots. This corresponds to 
a boost in the velocity of the corresponding pixels as shown in 
the in the distribution of [9*,,^'™'=] compared to [9*^^"™]. 

In the fourth row down in Figure 14 we find that even 
with the inclusion of the over dense region from the previ- 
ous test coupled with this additional under dense region, the 
final reconstruction shows a very similar result reconstructed 
velocity field where no additional observational effects have 
been embedded. With this very powerful test, we have shown 
that the Sy prior is once again not dominating the final re- 
constructed velocity field, and moreover, is indicating that the 
power spectrum (P^) prior and the residual signal in the ob- 
served velocity field V°^'^ is playing a key role in leading the 
reconstruction process to an unbiased result. 



So far, we shown the power of our method to successfully 
reconstruct the peculiar velocity field at redshifts out to z < 
0.3. In this section we show how the algorithm can be adapted 
to reconstruct the Newtonian gravitational potential 5*. 

From our definitions of the 5* and d'i', fields given in Equa- 
tions 32 and 34 respectively, one might suppose that trans- 
forming the final reconstructed velocity field into a 5* field is 
trivial. We can perform the integral by operating In Fourier 
space, 

nk) = ^ (43) 

However, there is a = surface which encodes a set of un- 
known radial bulk flow modes. We therefore need to take a 
different approach to the ^ reconstruction. 

7.1. Method 

The reconstruction algorithm for the ^ field is a modi- 
fied version of the velocity reconstruction method detailed in 
§ 5.3. We discuss the necessary modifications in relation to 
the chart in Figure 7 and the detail in § 5.3. 

Step 1: For our 5' reconstruction we set up the A^, number of 
initial 5'"]5';f' fields in the same way as we originally set up the 
trial velocity fields, so that for the n,/, initial field 



(44) 



where, ^^^fj is the true ^ field computed from the galaxy 
over density 6g as defined by Equations 3 1 and 32. The quan- 
tities and are gaussian random deviates with mean zero 
and variance of unity and / is the complex part. The parameter 
7] is set to the value used in the velocity reconstruction. 

Step 2: We now perturb the initial 5* fields in the same way 
as for the velocity reconstruction such that 

^TpAn) = *FFr(n > l) + (x" + iy")VP^(l^'\ix (45) 

where, as previously, x" and y" are uniform random deviates 
between [-0.5,0.5] and / is the complex pait. P*(A:)'™'' is the 
power spectrum of "i^pfj- For each n trial ^ reconstruction 
we can easily compute the coiTesponding velocity field by 



7. RECONSTRUCTING THE GRAVITATIONAL POTENTIAL 



y;?ry(n) = f5?^y(n)/k^ 



(46) 
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Figure 13. Artificially embedding a large structure into the ^P^"' field and allowing it to propagate into the 5v prior term. In both sets of panels, 
the top row shows the field as computed from the galaxy over-density 5g. The second row is the resulting gradient of that potential along 
the z-direction as denoted by d^z- The third row is the distribution between d'i'z and V""'^. The fourth row is the final reconstructed velocity 
field, V"^'^, for that particular slice in redshift. Finally, the bottom row is the distribution of V"''^ and V""*^. The red lines on the scatter plots 
indicate the least squares fit. For the left panel set at 0.135 < z < 0.141, the left-hand columns are the results where no artificial mass has 
been added to the simulations. In the right-hand column a gaussian mass has been added across 4 slices in z to the 'i' field, centred at X and Y 
pixel coordinates [11,12], The effect of this on the d^z field is then shown. In the right-hand panel set we can see this mass continuing in the 
subsequent redshift slice at 0.141 < z < 0.146. In both cases we can observe how this mass affects the resulting d'i'z field in the distributions 
of [a*-,!''™"] in the second row. 
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Figure 14. Similar test to that shown in Figure 13, where now we embed an under-density in the maps. 
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Step 3: For each trial reconstruction ^fpj{n) we compute 
its corresponding power spectrum Pl%(n). We also convert 

each Vffjin) to real space for evaluating in the log probability 
equation. 

Step 4: The log probability expression we use is now 
Xin) = }_^ 



Wpix 

/■=i 



l2 



a 



rpi 

1=1 



(«)■ 



^p*(k) 



E 



[y^'='=(n)-a^'(n)]2 



(47) 



Equation 46 shows that we must take the derivative of the trial 
reconstructed ^^ffj- field. This can lead to a noisy estimate of 
the resulting Vpfj given our choice of pixel scale. To min- 
imise this we have added an extra constraint in the form of 
the fourth term in Equation 47. That is, we retain the J,, prior 
from our velocity reconstruction method. This prior will com- 
plement the second term by guiding it toward less noisy solu- 
tions. Thus, for each iteration we evaluate the same velocity 
term as before shown as the first term in Equation 47. How- 
ever, we are now comparing the 'i' field and its corresponding 
power spectrum as shown by terms 2 and 3 respectively. 

Step 5: As before we determine the minimum log probability. 
If this corresponds to n = 1 we resample from the parent sam- 
ple in the next iteration of the code, and allow the pixel vari- 
ation parameter Ap,^ to reduce for the new set of trial fields. 
However, if the minimum log probability corresponds to value 
n 7^' 1, we set this to be the new template of trial velocity fields 
and return to Step 2. Steps 2 to 5 are repeated until conver- 
gence has been achieved. 

7.2. Results 

With our modified algorithm we now explore the ^ recon- 
struction for one realisation of the 0.08 < z < 0.18 cube simu- 
lation. For this test reconstruction we adopt in Equation 47 
as that given by Equation 27. The other input regularisation 
parameters are shown in Table 2. Figure 1 5 shows an exam- 
ple of spatial reconstruction compared to the true 5* field at 
the redshift slice 0.107 < z < 0.113. In the top panels we 
show the true gravitational potential field ^f*^"' for that slice. 
The second row shows the final reconstruction map, and 
the third row is the residual map e=(5'°^' Finally, the 

bottom left panel is the distribution of [ and the 

bottom right panel is the observed velocity distribution V°^^ 
compared to the true velocity distribution y'^'"'. In both of 
these panels we show the correlation r between the variables 
and the fitted slope a (also indicated by the red line). As with 
the velocity reconstructions we assume a ct = 0.3 error on the 



observed distances which propagates into our observed veloc- 
ity map. 

As we can see, the results are very encouraging. The recon- 
structed gravity map shows very good agreement with that of 
the true map with a strong correlation of r = O.823;I;0[J|4. 
As already observed in our velocity reconstructions, we find 
a bias in the slope of this distribution, here found to be 
a = 0.793 ± 0.029 which is also apparent in the residual map. 
Nevertheless, our algorithm reconstructs well the main struc- 
tures of the underlying true field. In Figure 1 6 we plot the cor- 
relation coefficients for each reconstructed slice as a function 
of redshift. We can see that out to z < 0.12 the reconstruc- 
tions show strong correlation with 4'°'''. Beyond this, there is 
an overall downward trend, reaching a minimum correlation 
of r = 0.460:1:0 at z = 0.155. While the overall main struc- 
tures are reconstructed, the amplitudes of the potentials seem 
to be offset compared to that of the true field; the reason for 
this will require further investigation. 

8. PROSPECTS WITH CURRENT SURVEY DATA 

Our simulations have allowed us to explore the potential 
of our reconstruction method as a cosmological probe that 
could be applied in future surveys, where the number density 
is sufficiently high. In this section we examine whether our 
approach is viable with the number density and distance esti- 
mator quality of current data. In recent years there have been 
a number of groups working toward building more compre- 
hensive catalogues of peculiar velocities. For example. Fun- 
damental Plane studies by Hyde and Bemardi (2009) have 
measured ^50,000 early type galaxies and more recently 
Magoulas et al. (2012) published 10,000 near-infrared early 
types in the 6dF Galaxy Survey out to z < 0.055. Current 
Tully-Fisher relation studies by Mocz et al. (2012) have mea- 
sured in excess of 25,000 Sloan Digital Sky Survey (SDSS) 
galaxies at 0.045 < z < 0.085. Recent work shows the emer- 
gence of another class of scaling relation to derive distances 
that relates mass and radius of early-type galaxies (see e.g. 
Chiosi et al. 2012; Bernardi et al. 2011), which may produce 
more precise velocity measurements than current DIs. 

Throughout this work our simulations contained ^300,000 
objects within a box that has the dimension of 384 h~' Mpc 
on the side which has then been shifted to be located at 
0.08 < z < 0.18 and 0.18 < z < 0.27. In Table 1 we sum- 
marised that this corresponded to respective number densities 
of ^ 140 and ^ 341 objects per square degree respectively. 
When we compare this to an approximate calculation of the 
SDSS Legacy data over the redshift range of 0.05 ^ ^ 0.15 
we find a number density of ^ 66 galaxies per square de- 
gree. Scaling this number to our 0.08 < z < 0.18 simulated 
cube would equate to approximately 144, 660 galaxies i.e. 
containing less that half the number of objects in our previ- 
ous analysis. Therefore, to test our algorithm on data with 
a more comparable signal-to-noise to what a current survey 
like SDSS could offer, we reduced the number of objects in 
each simulated cube. This was achieved by uniformly ran- 
domly removing a percentage of the parent data set to a level 
that may reflect a number density one might expect of com- 
bined measured distances estimated with both Tully-Fisher 
and the Fundamental Plane in SDSS. In both cubes we reduce 
the parent catalogues to a conservative number of ^ 70,000 
objects. This corresponds to a number density of ^ 32 and 
^ 79 galaxies per square degree for the 0.08 < z < 0.18 and 
0.18 < z < 0.27 respectively. We now apply our reconstruc- 
tion technique with the same parameter settings as before; the 
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Figure 15. Example of the final reconstruction of the ^' field for the 0.08 < z < 0.18 simulated galaxy cube. In this example we show the 
0.107 < z < 0.113 reconstructed slice. The top row panels show the true 4' field '^/^''K The second row shows the final reconstructed field 
vp''^'^ map and the third row is the residual map e between vp'^''' and *£'"^'^. Finally, the bottom-left panel is the distribution of [vl>'^'^'^,\l>^"'], and the 
bottom-right panel is the observed velocity distribution V"'"' compared to the true velocity distribution V'""^. In both of these panels we show 
the correlation r between the variables and the fitted slope a (also indicated by the red line). 
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Figure 16. The correlation coefficient values r for the final recon- 
struction of the VE" field for the 0.08 < z < 0.18 simulated galaxy 
cube. The top panels show the correlations between the residual 
maps e and '5°"' (shown in red) and 4'"^'^ (shown in blue). The bot- 
tom panel shows the correlation coefficient values between and 
the reconstruction 'i/"^'^ within each of the 16 slices in redshift of the 
reconstructed cube. 

results are summarised in Figures 17 and 18. 

The format of Figure 17 is the same as Figure 15. The 
left-hand panel set shows results for the 0.169 < z < 0.174 
slice and the right-hand panel set is for the 0.217 < z < 0.222 
slice. In both examples we see that, while the resulting ve- 
locity reconstruction is noisier than in previous tests, we can 
still achieve a good coiTelation and resolve the main struc- 
tures. Figure 18 shows the coiTelation coefficient values for 
all slices in this test with their associated residual coiTelations. 
For comparison we have superimposed the results from our 
earlier test in § 6 (Figure 10) shown as grey points. We can 
see that the result for this lower number density of objects fol- 
lows the same trend as our previous results in § 6. However, 
as expected there is a general systematic shift toward lower 
correlation values. Nevertheless, this is an encouraging result 
demonstrating that it may be possible to apply our technique 
existing data from e.g. SDSS, if we were to combine veloci- 
ties measured from both the Fundamental Plane and the Tully- 
Fisher relation. At the very least, forthcoming surveys should 
provide the required s/n level to reconstruct the gravitational 
potential from peculiar velocities. 

9. DISCUSSION & FUTURE WORK 

In this study we have attempted a first step towards a re- 
construction of the gravitational potential on cosmological 
scales, through the use of measured peculiar velocities as de- 
rived from cuiTent distance indicators. Previous methods have 
been centred around velocity reconstruction within the local 
Universe. Our motivation is to look toward the possibility of 
pushing to cosmological depths where such maps of the ve- 
locity and gravity field could be used as a cosmological probe 
ofGR. 

9.1. Velocity reconstruction 



For our velocity reconstruction algorithm we can sum- 
marise the approach and final results as follows: 

• We created a simulated observed velocity catalogue in- 
cluding distance indicator noise, drawn from a lognor- 
mal distribution motivated by the selection bias effect 
known as Inhomogenous Malmquist Bias. We have 
greatly simplified this effect by applying the errors ho- 
mogeneously and not allowing for density variations; 
we plan to address this in future versions of the method. 
We chose a conservative estimate for the statistical er- 
ror on the distances of cr^ = 0.3, which is greater than 
cuiTent error estimates from distance indicators. 

• Our simulated galaxy cube had a scale of 384 Mpc 
on the side and a total of 303, 575 objects. The cube 
was then pixelised to create a smoothed true velocity 
field y'™ and a velocity field that represents the noisy 
observed data, y"^*^. The cube was binned equally in x, 
y and redshift ZcosWith a transverse pixel scale of (24 
/i"' Mpc)^. For our analysis we consider the simulated 
galaxy cube at two redshift ranges: 0.08 < z < 0.18 and 
0.18 <z< 0.27. 

• Key to our method is the inclusion of two prior terms 
which guide the velocity reconstruction. The first prior 
term, featuring f^, regulates the reconstruction on var- 
ious scales by comparing the true velocity power spec- 
trum Pl'"''(k) to the one computed for the trial recon- 
struction Pl^'^(k). The second prior term, (5,„ compares 
the velocity field 9^- derived from the galaxy over- 
density 5g with the trial reconstructed velocity field 

yrec^ 

• The results showed that overall we can recover the ve- 
locity field to a high degree of correlation with y'™'=. 
For the cube at the 0.08 < z < 0.18 range, we found 
that we could achieve a consistently high correlation 
coefficient value of r[y'™'^,y'^'^'^] ^ 0.7 for a given 
slice in redshift with a high of r = 0.802 ± 0.012 at 
0.129 < z < 0.135. For the 0.18 < z < 0.27 cube we 
observed a similar trend in the final reconstruction with 
ahighof r = 0.808![i:j;j^ at 0.213 < z < 0.219. 

• We have also applied the algorithm to the same sim- 
ulation but with a much reduced number density of 
objects. By removing ^80% of the objects in the 
0.18 < z < 0.27 cube, we model a number den- 
sity that could possibly be achieved by an SDSS type 
survey if all sources of distance indicator information 
are utilised. From a total 306,896 objects in our par- 
ent catalogue, we were left with 70,488 objects within 
the same volume. Our results showed a similar trend to 
that of the original cube but with a reduced r[y'"'"^,y'""']. 
For example, for a slice at 0.174 < z < 0.180 we found 
r = 0.779 ±0.017 for the cube with 306,896 objects. 
However, with our reduced number density we found 
r = 0.606!j] for the same slice. Similarly, at 0.264 < 
z < 0.270 we find the constrained r = 0.752!°;g2i f^^^ 
to r = 0.588;!;Q yg2 for mock sample. 

9.2. Gravitational potential reconstruction 

Our method to reconstruct the gravitational potential re- 
quired a slight modification: 
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Figure 17. Resulting velocity reconstruction for SDSS type survey. In both panel sets we have reduced the number of galaxies in our simulated 
cubes from 300,000 to 70,000 to reflect a SDSS type survey. In the left-hand panel set we have selected the 0.169 < z < 0174 slice from 
the first simulated cube which ranges 0.008 < z < zO.18. The top row is the true velocity field V""*^, the second tow is the final reconstructed 
velocity field V^"^*^ and the third row is the residual map e between and V'"^'^. The bottom-left panel shows the [V"™'=,V"'^'^] distribution and 
the bottom-right panel is the [V'™'^,!/"''*] distribution. In both cases the correlation coefficient r is shown as well as the slope a also shown as the 
red line. The right-hand panel set shows the 0.217 < z< 0.222 reconstruction from the second simulated cube which ranges 0.18 < z < zO.27. 



• Instead of allowing the velocity field to evolve through 
the reconstruction process, we instead allow the ^I* field 
to be perturbed, converting then to a velocity field for 
evaluation in the likelihood. The new prior terms in 
and P^{k) are analogous to the 6y and P^^(k) pri- 
ors of the velocity reconstruction; with the S<s/ prior we 
construct a gravitational potential field computed di- 
rectly from the galaxy distribution via the Poisson equa- 
tion. We referred to this field as the true potential ^s'*' 
which we then compare with each trial reconstructed 
^j,iec rpjjg prior then compares the power spec- 

trum computed from the ^s^' field to the power spec- 
trum measure from ^f""^. In addition we retained the 
(5,. prior from the velocity reconstruction method as this 
allowed a further constraint to reduce the noise when 
computing the gradient of the potential d'i'-. 



• When applied to the 0.08 < z < 0.18 cube 
we found a fairly constant correlation coefficient of 
r[^r'<=<:,^'gai]=0.854 between 0.08 < z < 0.12, after 
which we saw the lowest r = 0.456;|^()9j occurring in 
the 0.152 < z < 0.157 slice. Beyond this, the recon- 
struction seemed to improve reaching r = 0.748^o in 
the final slice at 0.174 < z < 0.180. 

Overall, our results appear promising for peculiar velocity 
fields playing a more central role for future tests of gravity. 
Clearly, there are a number of aspects regarding this analysis 
that needs to be improved upon, and we will describe these 
below. 

9.3. Considerations for the future 

Code optimisation: There is much scope for optimisation 
of this algorithm, so future investigations should study 
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Figure 18. Correlation coefficient values for the final reconstruction 
of the peculiar velocity field for a SDSS type survey sample. The 
bottom panels show the correlation between r[l/°''\V'™'^] (red points) 
for each redshift slice in our simulated galaxy cubes; 0.09 < z < 0. 1 8 
(left panels) and 0.18 < z < 0.270 (right panels). The black points 
show the correlation of the true velocity field i.e. r[V'^'^'^,V""'^] and the 
grey points show r[V'^'^^,V'""^] from Figure 10. The top panels show 
the respective residual correlations between the V""'' (black points) 
and V'"" (blue points), where e is the residual, e = V'™" - V""'. 

the effects on the reconstruction of with varying sizes of 
smoothing scales together with more sophisticated smoothing 
kernels applied to the input data. It will also be valuable 
to optimise the way in which we model our regularisation 
parameters for the priors. In the current version they are fixed 
quantities and so not allowed to adapt with each iteration of 
the algorithm. 

Mock catalogues and modelling selection biases: Here 
we have applied a simple approach to modelling the effect 
of Inhomogeneous Malmquist bias, by drawing distance 
errors from a lognormal distribution and applying it to 
our data homogeneously. Therefore we have not taken 
into account the density variations along the line of sight 
which are responsible for this type of bias. Clearly, when 
adapting our algorithm to real survey data this will have to be 
investigated more carefully. Whilst our estimates for a SDSS 
type survey are crude; we intend to create more sophisticated 
mock galaxy catalogues which will encapsulate all typical 
selection effects for specific surveys. Nevertheless, this work 
represents a good first attempt towards our goal of being able 
to apply this reconstruction approach to real survey data. 

Kaiser Effect & Fingers of God: In the cuiTent analysis 
we have neglected possible adverse affects on the recon- 
struction that may arise from redshift space distortions. On 
large-scales these take the form of the Kaiser effect where 
objects fall coherently towards over dense regions and away 
from under dense regions; this phenomenon is the physical 
basis for our reconstruction. But on small scales, non-linear 
effects that arise in redshift space cause the so-called Fingers 
of God effect where we observe a smearing of objects along 
the line of sight. A full treatment including Fingers of God 



will be a further improvement reserved for future work. 

Cosmic Variance: We are also currently testing our algorithms 
on a Hubble-like Volume simulations (see Warren et al. 2012) 
in an attempt to quantify any potential uncertainties due to 
cosmic variance (e.g., Somerville et al. 2004). 
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